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ABSTRACT 

We observed a secondary eclipse of WASP-33 b quasi-simultaneously in the optical (~0.55 pm) and the near-infrared (~1.05 pm) 
using the 2x8.4 m Large Binocular Telescope. WASP-33 is a d Scuti star pulsating with periods comparable to the eclipse 
duration, making the determination of the eclipse depth challenging. We use previously determined oscillation frequencies to 
model and remove the pulsation signal from the light curves, isolating the secondary eclipse. The determined eclipse depth is 
AF = 1.03±0.34 parts per thousand, corresponding to a brightness temperature of Tb = 3398 ± 302 K. Combining previously pub¬ 
lished data with our new measurement we find the equilibrium temperature of WASP-33 b to be Tb = 3358 ± 165 K. We compare all 
existing eclipse data to a blackbody spectrum, to a carbon-rich non-inverted model and to a solar composition model with an inverted 
temperature structure. We find that current available data on WASP-33 b’s atmosphere can be best represented by a simple blackbody 
emission, without the need for more sophisticated atmospheric models with temperature inversions. Although our data cannot rule out 
models with or without a temperature inversion, they do confirm a high brightness temperature for the planet at short wavelengths. 
WASP-33 b is one of the hottest exoplanets known till date, and its equilibrium temperature is consistent with rapid reradiation of the 
absorbed stellar light and a low albedo. 
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1. Introduction 


Hot Jupiters offer us a great astrophysical laboratory for atmo¬ 
spheric science under extreme conditions as their atmospheres 
reach temperatures above 1000 K. Secondary eclipse measure¬ 
ments easily separate the stellar from the planetary lights and 
are therefore usually used to probe exoplanet atmospheres. Such 
measurements are most straight forward in the near infrared 
(NIR), where the planet-to-star contrast ratio is large (see e.g., 
Deming et al. 2007 [ |Christiansen et ah |2010[ Knutson et al. 


201 l|l. The first detection of secondary eclipses were made from 


the space using the infrared Spitzer Space Telescope (Charbon- 
neau et ari|2005[ Deming et al.||2005| l. Nonetheless, the plane 
tary atmospheric composition and temperature structure has also 
been detected from the ground (see e.g., Snellen et al. 2008[ 
|Bean et al.|2010| l. Early studies suggest the existence of a tem¬ 
perature inversion layer caused by strong stellar absorption of 
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certain molecules high up in the exoplanet atmospheres. Eor hot 
Jupiters, theoretical predictions place TiO and VO as favorite 
candidates (see e.g., Hubeny et al.||2003|l. However, this is still 


under debate ( [Spiegel et al. 2009| l. |Eortney et al.| ( |2008| l proposed 
a differentiation of hot Jupiters into two populations divided by 
the amount of stellar irradiation. Eor these two families, the tem¬ 
perature profiles would be dominated by the state of TiO and VO. 
However, more recent observations have not been able to clearly 
show the existence of these two states ( [Knutson et al.|[2()10t 
[Crossfield et al.[2012| l. A good example is WASP-12b, for which 
no inversion has been yet detected ( [Crossfield et al.[2012| l despite 
being one of the most irradiated and hottest exoplanets known 
(Teq ~3000 K, Hebb et al. 2009 Lopez-Morales et al. 2010 l 
A s observing techniques and data reduction strategies improve 
with time, past claimed detection of temperature inversions are 
being overruled. Eor example, HD 209458b was observed us¬ 
ing Spitzer Space Telescope by [Knutson et al.[ ( |2008[ l. The four 
bandpass observations were represented by a thermal inversion 
in the exoplanet atmosphere (see e.g., [Line et al.[2014) l. However, 
[Diamond-Lowe et al.[ ( [2014| ) used state of the art techniques and 
models to perform a complete and self-consistent analysis of the 
data, finding no temperature inversion after all. To date, there 
doesn’t seem to be a clear determinant for which planets should 
present (or not) temperature inversions. Measuring the structure 
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of one of the hottest exoplanet atmospheres could help to clarify 
the situation. 

This work focuses on the transiting hot Jupiter WASP-33 b 
( [Collier Cameron et al.||2010| l. With a wavelength-dependent 
brightness temperature of 3620 K, WASP-33 b is one of the 
hottest exoplanets known to date ( Smith et al.|201l| l. WASP-33 
(HD 15082) is a (5 Scuti star, osci lating with amplitudes in the 
milli-magnitude regime and periods that are comparable to the 
transit duration ( jHerrero et al.|20lT] l. WASP-33 b orbits the star 
each ~1.22 days in a retrograde orbit ( [Collier Cameron et al.j 
2010[ l. Showing an unusually large radius, WASP-33 b belongs 
to the class of anomalously inflated exoplanets ([Collier Cameron[ 
jet al.|2010j l. 

The first secondary eclipse was detected by Smith et al. 


( [201 l| l at 0.91 pm. [Deming et al.[ ( |2012[ ) observed two events 
in the Ks, band (2.15 pm), and one secondary eclipse each at 
3.6 pm and 4.5 pm using Warm Spitzer. The authors found the 
data to be consistent with two different atmospheric models, 
the first one with an inverted temperature structure and a solar 
composition, and the second one with a non-inverted tempera¬ 
ture structure and a carbon-rich composition. Furthermore, [de| 
Mooij et al. (2013 1 re-observed WASP-33 b’s eclipses in the Ks, 
finding consistency with Deming et al. ( 2012| l and determined 
WASP-33 b’s effective temperature to be around 3300 K. Re¬ 
cently, Haynes et al. ( [2015| l observed WASP-33 during two oc- 
cultations using the Hubble Space Telescope between 1.13 and 
1.63 pm. The authors represented their observations by a day- 
side atmosphere with a temperature inversion and an oxygen- 
rich composition with a slightly sub-solar abundance of H 2 O. 
Excluding this work, previous studies on WASP-33 b’s emission 
were observed in a broad wavelength range. At that time the au¬ 
thors didn’t count with a data-independent pulsation model nor 
with the characterization of the pulsation spectrum of the host 
star. In von Essen et al. ( 2014[ l we presented a two and a half 
years long follow-up campaign of WASP-33. We collected 457 
hours of in- and out- of transit data that were used to characterize 
the pulsation spectrum of WASP-33 with unprecedented preci¬ 
sion, finding 8 statistically significant frequencies with periods 
comparable to the transit duration (~2.7 hours). In this work, we 
observed WASP-33 during secondary eclipse at 1.05 pm using 
the Large Binocular Telescope. We make use of the identified 
pulsations to “clean” the photometry and determine the eclipse 
depth of WASP-33 b’s secondary eclipse. In this work. Section]^ 
presents the observations and the data reduction, Section|^a de¬ 
tailed description of the model parameters and fitting procedures. 
Section]^ the secondary eclipse measurements in the context of 
the atmospheric structure of WASP-33, and Section]^ our final 
remarks. 


2. Observations and Data Reduction 

We observed WASP-33 on October 12/13, 2012, during a sec¬ 
ondary eclipse with the 2x8.4 m Large Binocular Telescop^ 
(LBT) in binocular mode. In particular, we employed the NIR 
spectrograph LUCI ( [Seife rt et al.[[201 0|) simultaneousl y to the 
Large Binocular Cainera^(LBC Giallongo et al. 2008 1 . Unfor¬ 
tunately, due to technical problems the spectroscopic LUCI time- 
series lacks the pre-ingress part of the light curve. Therefore, our 
work focuses on the photometric data only. 

The LBC consists of two wide-held cameras mounted on the 
prime focus swing arms of the LBT. Photometric data were ac- 


* www.lbto.org 
^ lbc.mporzio.astro.it 


quired using the LBC Red camera, which is optimized for ob¬ 
servations between 0.55 and 1 pm. The camera has four charge 
coupled devices (CCDs), each one of them with its own gain 
and read-out noise characteristic 0 To minimize read-out time 
we used sub-frames of the CCDs with the brightest stars within 
the held of view centered on each chip. The size of the hnal 
sub-frames were approximately 1 x3 arcmin, which reduced the 
read-out time from 30 to 15 seconds. 

During the night the observing conditions were stable and 
the sky was clear. Seeing varied between 0.4 and 1.1 arcsec. 
The photometric data were acquired using the Johnson-Cousins 
V hlter centered at 0.55 pm and a Y hlter centered at 1.05 pm, 
the latter with a full-width at half-maximum of 0.2 pm. As 
WASP-33 is relatively bright (V~8.3), we defocused the tele¬ 
scope to avoid saturation. The transmission function of the Y- 
band peaks around 95%, but the quantum efficiency of the CCD 
drops down to 10% in the hlters wavelength range. This was also 
useful to avoid saturation. WASP-33 is located in a sparse stel¬ 
lar field. Therefore, the defocusing technique did not produce 
overlapping of the stellar point spread functions, except for the 
optical companion identified by Moya et ar| ( |201 l| l. A discussion 
on third-light contribution is addressed in Section [33| Exposure 
times were of 0.74 seconds for the V-band and 3.23 seconds for 
the T-band. Every 3 images the filters were switched between the 
two bands. Since the readout time is shorter than the time that the 
filter wheel requires to go from one filter to the other (~20 sec), 
the sampling between consecutive images is each ~18 seconds, 
while the sampling between consecutive groups of 3 images of 
the same filter is each ~2.3 minutes. 

We used the IRAE’s package ccdproc for bias subtraction 
and the flat fielding, and the IRAE task apphot for aperture pho¬ 
tometry. We measured fluxes of WASP-33 and three reference 
stars, namely BD-h36 493, BDh- 36 487, and BDh- 36 488, us¬ 
ing various aperture radii. We then fixed the aperture radius so 
that the scatter of the differential light curve is minimized. The 
chosen aperture radius was of 74 pixels, equivalently to 16 arc- 
sec considering the plate scale of the LBC Red (0.2255 "/pixel). 
Next we divided WASP-33’s measured flux by the average of the 
three reference stars producing the final differential light curves 
(see Pigure[2and Eigure]^. Since these are influenced by stellar 
pulsations, we estimated the scatter averaging the standard de¬ 
viations calculated within 15 minutes bins, where the deforma¬ 
tions caused by WASP-33’s intrinsic variability can be neglected. 
The derived uncertainties are ~0.8 parts-per-thousand (ppt) for 
the V-band and 0.6 ppt in the T-band, four and three times the 
expected photometric precision assuming only photon noise (0.2 
ppt), respectively. We believe that the difference between the 
achieved photometric precision and the theoretical photon-noise 
precision is caused by scintillation noise due to the short expo¬ 
sure times. To ensure that the photometry is based on proper ref¬ 
erences, the three reference stars were plotted against each other 
to check for constancy. 

During the observations, the LUCI (spectroscopic) side of 
the LBT was the one guiding the telescope. The LBC side fol¬ 
lowed, however, its pointing stability was rather poor. After 3.5 
hours of observations WASP-33 and the reference stars drifted 
~200 pixels in both x and y directions close to the edge of the 
read-out windows and flux was lost. These ~0.5 hour of ob¬ 
servations were omitted from our analysis. WASP-33 b eclipse 
lasts 2.7 hours. Therefore, our observations account for ~25 
minutes of off-eclipse data before eclipse begins and ends, re¬ 
spectively, and full eclipse coverage. As a final step the time 
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stamps were converted from Julian Dates to Barycentric Julian 
Dates (BJDtdb) using the tools made available by|Eastman et al. 
( |20T01 i. 


3. Results 

3.1. Analysis of the optical light curve 

The data presented here comprise quasi-simultaneous observa¬ 
tions during secondary eclipse of WASP-33 b around the V and 
Y bands. The predicted planet-star flux ratio in the T-band is 
0.2 ppt, four times lower than the accuracy of our measurements. 
Therefore, we can neglect the planet imprint and use this band 
to measure the stellar pulsations, and most specihcally to tune 

Ivon Essen et al.|2014|l. 


their current phases (see Phase Shifts in 
Particularly, our model for the light contribution of the stellar 
pulsations consists of eight sinusoidal pulsation frequencies with 
corresponding amplitudes and phases. Hence, to reduce the num¬ 
ber of 24 free parameters and the values they can take, we use 
prior knowledge about the pulsation spectmm of the star that 
was acquired durin g |von Essen et al. ( 2014| l. As the frequency 
resolution is 1/AT ( |Kurtz|1983| l, 3.5 hours of data are not suffi¬ 
cient to determine the pulsations frequencies. Therefore, during 
our fitting procedure we use the frequencies determined in |von| 
jEssen et al.| ( |2014| l as starting values plus their derived errors 
as Gaussian priors. As pointed out in von Essen et al. ( 2014| l, 
we found clear evidences of pulsation phase variability with a 
maximum change of 2x10“^ c/d. In other words, as an exam¬ 
ple after one year time a phase-constant model would appear to 
have the correct shape with respect to the pulsation pattern of the 
star, but shifted several minutes in time. To account for this, the 


eight phases were considered as fitting parameters. The von Es- 
|sen et al.] ( |2014| l photometric follow-up started in August, 2010, 
and ended in October, 2012, coinciding with these LBT data. 
We then used the phases determined in von Essen et al. ( |2014| l 
during our last observing season as starting values, and we re¬ 
stricted them to the limiting cases derived in Section 3.5 of|von| 


the blue wavelength range. Therefore the amplitudes estimated 
in [von Essen et al.|P014| l, listed in Table[T] are used in this work 
as fixed values. This approach would be incorrect if the pulsa¬ 
tion amplitudes would be significantly variable (see e.g., |Bregei^ 
|et al.|2005] l. Nonetheless, the short time span of LBT data, and 
the achieved photometric precision compared to the intrinsically 
low values of WASP-33’s amplitudes, make the detection of any 
amplitude variability impossible. 

Eor the fitting procedure we used the photometric errors that 
IRAF provides, but enlarged to meet the standard deviation of 
the data. Our model to the V-band light curve requires three 
additional free parameters to account for airmass and color- 
dependent extinction. Since the pulsations impede any visual 
inspection of how these effects impact on the photometry, the 
number of parameters was determined through minimiza¬ 
tion. To derive the 16-1-3 parameters for our model we used a 
Markov-Chain Monte Carlo (MCMC) approach. Our MCMC 
calculations make extensive use of PyAstronomjj^ a collection 
of Python routines providing a convenient interface for fitting 
and sampling algorithms implemented in the PyMC (|Patil et 


2010i and SciPy (Jones et al. 200 1| packages. After 1x10^ 
MCMC iterations and a burn-in of the initial 2x10^ samples, we 
computed the mean and standard deviation (1-cr) of the posterior 
distributions of the parameters and used them as best-fit values 
and uncertainties, respectively. Our results for the frequencies 
and phases are summarized in the sixth and seventh columns of 
Table The small errors on the fitted frequencies are caused due 
to the priors imposed and due to the reduction of the fitting pa¬ 


rameters compared to the ones in von Essen et al. ( 2014) l. In this 
work the eight amplitudes are considered as fixed parameters, 
making as consequence the statistical errors to appear smaller. To 
test if the MCMC chains were well-mixed and converged we fol¬ 
lowed two approaches. Firstly, we started the MCMC chains us¬ 
ing different initial values. Secondly, we divided complete chains 
into sub-chains (normally dividing the 1x10® iterations into 4) 
confirming afterward that the best-fitting values obtained from 
the posterior distributions of each sub-chains were consistent to 
each other. Finally, we confirmed convergence by visually in¬ 
specting the traces and histograms of all the fitting parameters. 
Figure [T] shows the LBT optical data over-plotted to the best-fit 
model (top) along with the residual light curve (bottom). 

Our analysis over the V-band data was concluded trying to 
estimate an upper limit for the eclipse depth. For this end, we 
added to our model budget a secondary eclipse shape (see Sec- 
tion |3.2.2 for details on secondary eclipse modeling). After sim¬ 
ilar iterations and sample burn-in, we found the best-fit eclipse 
depth to be AF = 0.26 + 0.48 ppt. computations between the 
two models favor the case where the secondary eclipse is not 
considered. As previously assumed, the derived uncertainty is 
too large to compute a reasonable 3-cr upper limit. 


Orbital phase 

- 0.06 - 0.04 - 0.02 0 0.02 


0.04 0.06 


Essen et al. (20141, rather than allowing them to take arbitrary 
values. The pulsation amplitudes in 6 Scuti stars are expected 
to be wavelength-dependent (see e.g. |Daszyriska-Daszkiewicz| 
2008 |l. Our follow-up campaign and these data were acquired in 


CD 

DC 



0.996 - 


- 2.0 - 1.5 - 1.0 - 0.5 0.0 0.5 1.0 1.5 2.0 

Hours from transit 

Fig. 1: LBT V optical data during secondary eclipse of 
WASP-33 b as a function of the orbital phase and the hours rela¬ 
tive to mid-eclipse. Top: Normalized flux. The continuous black 
line shows the best-fitting model (24-1-3 parameters). The error 
bars were enlarged to meet the scatter of the data. Bottom: resid¬ 
ual light curve. Vertical dashed lines indicate first to fourth con- 


http://www.hs.uni-hamburg.de/DE/Ins/Per/Czesla/ 
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tacts assuming a circular orbit (Deming et al. 20121. 
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Table 1: Comparison of derived pulsations between previous and current work. From left to right: Frequency number (PN), fre¬ 
quencies in c/d, amplitudes in ppt, phases in degrees and range of variability for the phases determined in |von Essen et al.1 ( |2014| l 
and used here as model parameters. Additionally, derived frequencies and phases {V,Y) found in this work and phase differences in 
degrees. Errors are in all cases at 1-cr level. 


PN 

V 

(c/d) 

A 

(ppt) 

(/v 

(°) 

A(/v 

(°) 

VV,flt 

(c/d) 

i/v,fit 

(°) 

VY,flt 

(c/d) 

0Y,fit 

(°) 

i/y-v 

(°) 

fi 

20.1621 ±0.0023 

0.95 ± 0.03 

241 ± 1 

[170:246] 

20.1636 ±0.0007 

243 ± 8 

20.1620 ±0.0012 

258 ± 12 

15 ± 15 

f2 

21.0606 ±0.0023 

0.93 ± 0.03 

122 ± 1 

[21:185] 

21.0616 ± 0.0007 

171 ± 7 

21.0577 ± 0.0089 

157 ± 13 

-14 ± 15 

f3 

9.8436 ± 0.0023 

0.79 ± 0.03 

195 ± 1 

[82:258] 

9.8445 ± 0.0009 

249 ±7 

9.8428 ± 0.0016 

274 ± 11 

25 ± 13 

f4 

24.8835 ± 0.0017 

0.42 ± 0.03 

115 ± 2 

[54:126] 

24.8835 ± 0.0004 

118 ± 7 

24.8835 ± 0.0009 

115 ± 11 

-3 ± 13 

fs 

20.5353 ± 0.0013 

0.71 ± 0.03 

212 ± 1 

[115:275] 

20.5329 ± 0.0010 

256 ± 10 

20.5356 ±0.0011 

274 ± 12 

18 ± 15 

f6 

34.1252 ±0.0027 

0.49 ± 0.03 

218 ± 2 

[152:230] 

34.1227 ±0.0004 

223 ± 7 

34.1268 ±0.0009 

230 ± 8 

7 ± 11 

f7 

8.3084 ± 0.0025 

0.63 ± 0.03 

194 ±2 

[42:204] 

8.3107 ±0.0010 

188 ± 10 

8.3062 ± 0.0012 

202 ± 16 

14 ± 19 

fs 

10.8249 ± 0.0030 

0.64 ± 0.03 

11 ±2 

[10:196] 

10.8236 ± 0.0010 

110± 8 

10.8256 ± 0.0061 

103 ± 13 

-7± 15 


3.1.1. Residual light curve analysis 

Pulsations not accounted for in our model would manifest them¬ 
selves as correlated noise. To quantify to which extent the resid¬ 
ual light curve is affected by it, we followed a similar approach 
as Gillon et al.| ( |2006[ ); |Winn et ah] ( |2008[ ), and [Carter & Winn| 


(2009|l. We started dividing the residual light curve into M bins 


of equal duration. Since our data are not equally spaced, we cal¬ 
culated a mean value N of data points per bin. If the data are 
affected by correlated noise, the sample standard deviation of 
the binned data, ctn, would differ by a factor /3 n from its the¬ 
oretical expectation (see e.g., von Essen et al.||2013| for an ex¬ 
tended mathematical description). For data sets free of corre¬ 
lated noise, /^n = 1 is expected. ^Sn’s smaller than 1 might occur 
due to statistical fluctuations and were neglected in our analy¬ 
sis. We estimated the amount of correlated noise, ff, averaging 
jSn’s that were calculated within bins of 0.5, 0.75, 1, 1.25 and 
1.5 times the ingress time. Its derived value, (3 - 1.05, allows us 
to neglect the contribution of correlated noise in our future anal¬ 
ysis. Finally, we searched for any residual pulsations computing 
a Lomb-Scargle periodogram ( |Lomb|1976[ Scargle 1 1982 1 [ZechT] 
jmeister & Kurster|2009| l from LBT-V’s residual data and found 
no signihcant peaks. Thus, our pulsation model satisfactory rep¬ 
resents the intrinsic variability of WASP-33. 


3.2. Analysis of the NIR light curve 

3.2.1. Wavelength-dependency of fitting parameters for 6 
ScutI stars 

The pulsations of 6 Scuti stars are known to be wavelength- 
dependent (see e. g., [Balona & Eversj |1999[ jPaszyhskaT] 
Daszkiewicz et al. |20()3| l. Thus, subtracting the F-model to 
the Y data would not properly account for WASP-33’s intrin¬ 
sic variability. In [von Essen et al.j ( |2014| l we tried to charac¬ 
terize the nature of the modes analyzing simultaneous multi¬ 
color photometry (see Section 3.6.3 of [von Essen et aL]p014| l 
for further details). However, we found no conclusive results. 
Consequently, an estimation of the V-to-F change of the am¬ 
plitudes and phases of the eight frequencies can only be ob¬ 
tained via theoretical modeling. To quantify these changes we 
used the Erequency Analysis and Mode Identihcation for As- 
teroSeismology code (EAMIAS), a collection of tools for the 
analysis of photometric and spectroscopic time series data (jZimaj 


2008 j l. To begin with, EAMIAS requires basic stellar parameters. 
Therefore, we adopted T^ff = 7430 + 100 K, log(g) = 4.3 + 0.2, 
and [Ee/H] = 0.1 + 0.2 ( [Collier Cameron et al.|201^ . The non- 
adiabatic observables were calculated using the code MAD 
( [Dupret[|2001[ [Grigahcene et al.|[2005[ ), while as model atmo- 
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spheres we used [Kurucz[ (19971’s grids. Eurthermore, we as¬ 
sumed the degree, (, and the azimuthal number, m, of the pul¬ 
sation modes to be low {m, f < 2, see e.g. [Aerts et'ar][2010| 
geometric cancellation). Erom EAMIAS we computed the am¬ 
plitude ratios and phase differences for each one of the eight fre¬ 
quencies relative to the V-band. Unfortunately, the wavelength 
coverage of EAMIAS goes up to 0.8 pm, i.e. not covering the Y 
wavelength range. Therefore, our estimates are an extrapolation 
of EAMIAS results, all of them summarized in Table As an 
example, Eigure shows EAMIAS results for vi - 20.1621 c/d, 
along with their respective extrapolations. Analyzing these the¬ 
oretical expectations some relevant conclusions can be drawn. 
Eirst, as it is expected for d Scuti stars the amplitudes of the 
pulsations are larger where the star emits most of its flux (see 
e.g., |Daszynska-Daszkiewicz|2008| l. Therefore, for WASP-33 it 
is expected that the amplitudes will decrease for redder wave¬ 
lengths. Second, the computed amplitude ratios for different { 
degrees would make the proper identification of the modes im¬ 
possible to achieve by means of these multi-color data. Indeed, 
for the eight frequencies that comprise our model the amplitudes 
are expected to scale down around 40-80% in the Y with respect 
to their V values. Considering the relative short duration of these 
observations, the photometric precision of the data, and the rela¬ 
tion between this precision and the amplitude of the pulsations, 
for low Avalues it is impossible to determine their geometry (see 
Eigure]^. Einally, the phases in V and Y are expected to be sim¬ 
ilar. All these aspects will be considered in our NIR pulsation 
model. 


Table 2: Maximum theoretically-expected changes for am¬ 
plitude ratios {Ay jAy min-Ay jAy max) phase differences 

{(pvxmin,rpv,Y,ma.x) between the V and F-bands for the eight pul¬ 
sation frequencies considering f <2. 


PN 

Frequency (c/d) 

^■^max 

A4>max C°) 

f'l 

20.1621 

0.3-0.7 

(-30,20) 

f2 

21.0606 

0.6-0.9 

(-20,25) 

fa 

9.8436 

0.5-0.8 

(-23,20) 

U 

24.8835 

0.5-0.8 

(-20,20) 

fs 

20.5353 

0.4-0.75 

(-25,20) 

f6 

34.1252 

0.5-0.7 

(-20,20) 

fv 

8.3084 

0.45-0.75 

(-30,20) 

fg 

10.8249 

0.35-0.75 

(-25,20) 
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1 = 0 - 1 = 1 - 1 = 2 



0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 

Wavelength (|xm) 


Fig. 2: Amplitude ratios (top) and phase differences bottom 
obtained with FAMIAS for the vi - 20.1621 c/d pulsation fre¬ 
quency, considering £ <2. Gray contours indicate the extrapo¬ 
lated area extended to the Y filter. Both V and Y filter transmis¬ 
sion functions are plotted in gray continuous lines. 


3.2.2. Planetary size versus secondary eclipse model 

As secondary eclipse model we used a simple step function. Its 
shape was limited by the first to fourth contacts ( |Sackett|1999| l. 
The latter were estimated using the semi-major axis a/R^, the 
orbital inclination i, the planet-to-star radius ratio Rp/Rs, the or¬ 
bital period P, and the mid-eclipse time Tq that were computed 
in |von Essen et al!] ( |2014[ ) accounting for the stellar pulsations. 
To was used to place the step function relative to mid-eclipse 
time, assuming a circular orbit ( |Deming et al.|201^ . The chosen 
values are summarized in Table[3] 


Table 3: Orbital and physical parameters of WASP-33 b com¬ 
puted in von Essen et al. (2014 1 . 


Orbital parameters Best-fit values and 1 

-cr errors 


von Essen et al. 

2014)) 

P (days) 

1.2198675 + 11x10 ° 

To (BJDtbd) 2455507.5222 + 0.0003 

a/Rs 

3.68 + 0.03 

i (deg) 

87.90 + 0.93 

Rp/Rs 

0.1046 + 0.0006 


Any potential change in the planetary radius as function of 
observing wavelength mainly manifests itself in a change of the 
eclipse shape. In addition, the times of contact might change by 
a small amount. We tested if our data is sensitive to this second 
order effect and if the use of one set of contact points for all 
wavelengths is justified. We simulated a secondary eclipse of a 
transiting planet 20% larger than the value measured in the opti¬ 
cal (Table[^, adding to it the time stamps and noise properties of 
our data. We find that a change in the eclipse duration does not 
improve the fit to the simulated data compared to a fit where the 
radius obtained in the optical is used. Our data quality is such 
that we do not need to take this second order effect into account 
and we proceed using the times of contact as determined in the 
optical. 


3.3. Fitting parameters and resuits 

We now use our knowledge on the pulsations gained from the 
analysis of the optical light curve in the analysis of the NIR light 
curve. Our model consists of an eclipse light curve described 
by the parameters from Table For all but the planet to star flux 
ratio AF, which we want to determine, we use the parameters and 
uncertainties given in the table as prior information. To describe 
the stellar pulsations we use the same model as in the optical 
plus the prior information on the eight frequencies, amplitudes 
and phases. All amplitudes are damped by a common factor AA 
relative to the optical. All parameters and their 1-cr uncertainties 
are given in Table [T] In particular, we find an eclipse depth of 
AF = 1.03 + 0.34 ppt in the T-band. 

The derived frequencies, phases, and the computed phase 
differences between both bands, are summarized in the last three 
columns of Table [1] Our results show that both cpy - (py and 
AA = 0.73 + 0.08 are consistent with their theoretical expecta¬ 
tions. Figure shows in the top panel the T-band light curve, 
in the middle panel the pulsation model removed isolating the 
planetary eclipse and in the bottom panel the residuals model 
data. On top of our best-fit secondary eclipse model we show 
|Deming et al.| ( |20T^ ’s predicted eclipse depths integrated in the 
T-band for the atmosphere of WASP-33 b when a carbon-rich 
non-inverted model (green) and a solar composition model with 
an inverted temperature structure (pink) are considered. 

We calculated the Pearson, P, and Spearman, S , correlation 
coefficients between AA and AF. Both show a low positive cor¬ 
relation of ~0.4. No such correlation was found between the 
eclipse depth and the phases nor frequencies. 


3.4. Robustness of our fitting approach 

Measured eclipse depths and uncertainties can depend sensi¬ 
tively on the choice of model used to fit the light curve. To con¬ 
firm that our final result is robust, we performed a number of 
additional fits. We compared the models using the Bayesian in¬ 
formation criterion, BIC = ;t' -H k In N, which penalizes the num¬ 
ber k of model parameters given N data points. The additional 
models for the stellar pulsations included the following fitting 
parameters: 


2 . 


Frequencies fixed to the best-fit values found in von Essen 
|et al.| ( |2OT4) l, scaling factor for the amplitudes, eight phases 
with Gaussian priors, secondary eclipse and detrending func¬ 
tion. 

Eight frequencies with Gaussian priors, eight amplitudes 
with uniform priors, eight phases with Gaussian priors, sec¬ 
ondary eclipse and detrending function. 


3. Frequencies fixed to the best-fit values found in von Essen 


et al. (2014 1 , eight amplitudes with uniform priors, eight 


phases with Gaussian priors, secondary eclipse and detrend¬ 
ing function. 


The resulting BIC values, along with the fitting parameters 
for each model and the derived best-fit eclipse depth for each 
model are summarized in Table As the table reveals, the BIC 
favors our approach. All eclipse depths are consistent within 1-cr 
errors. 


3.5. Potentiai third-iight contributions 


Moya et al. ( 201 l| l detected an optical companion to WASP-33 
at an angular separation of ~2”. This companion is included in 
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Fig. 3: Equivalently to Figure but around 1.05 pm. From top 
to bottom, WASP-33’s pulsations and secondary eclipse, the sec¬ 
ondary eclipse once the pulsations have been removed, and the 
residual light curve once the secondary eclipse is removed. Black 
continuous lines show best-fit models for the pulsations and the 
secondary eclipse. Pink and green eclipse models correspond to 


Deming et al. (2012 1 ’s predicted eclipse depths for inverted and 


non-mverted atmospheric structures integrated in the F-band. 


Table 4: Results for the four different pulsation models. From 
left to right, stellar pulsation model (Model 0 corresponds to the 
one carried out originally, while the remaining numeration fol¬ 
lows the previously detailed stellar models), the number of fitting 
parameters, k, the computed value between the observations 
and the best fitting models, the Bayesian Information Criterion 
(BIC) and the derived eclipse depth (ED). 


Model 

k 


BIC 

ED (ppt) 

0 

21 

25.3 

141.1 

1.03 + 0.34 

1 

13 

243.4 

315.1 

1.23 + 0.24 

2 

28 

27.3 

181.7 

0.98 + 0.45 

3 

20 

249.1 

359.4 

1.24 + 0.25 


our aperture and therefore reduces the measured eclipse depth 
relative to the undiluted case. Based on the color information 
obtained in the Jc, He, Kc, and Fll filters, [Moya et al.| ( |201 l) l es¬ 
timated that the probable physical companion of WASP-33 has 
an effective temperature of Teff = 3050 + 250 K. To compute the 
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contribution of this third light we used the approach detailed in 
Tingley p004)l. First, we represented the emission of the third 


body, F 3 , relative to WASP-33’s emission, Fw 33 , with the ratio 
of two black body curves F 3 /FW 33 of 3050 and 7430 K, respec¬ 
tively. Then, we converted the Amag’s of the Jc, He, Kc, and 
FII bands to fluxes and adjusted F 3 /FW 33 to Moya et al. ( |201 l[ )’s 
measurements with a diluting (scaling) factor, C, that accounts 
for the relative distances between the two stellar objects. Extrap¬ 
olating this to the T-band we found F 3 to be 714 times dimmer 
than WASP-33. Given our data we can ignore such a contribu¬ 
tion. Although strictly speaking stars are not blackbodies, the 
blackbody approximation was carried out in this case mainly 
due to two reasons. On one hand, WASP-33 is an A-type star 
and, therefore, is expected to have no strong absorption lines in 
the wavelength range of our interest. On the other hand, there 
is no knowledge about the metallicity nor the surface gravity of 
WASP-33’s companion, not allowing us to use appropriate stel¬ 
lar models without further speculation. Nonetheless, we caution 
the reader to carry out this approximation when a more com¬ 
plete knowledge of the spectral properties of the involved stars 
is known. 


4. Thermal emission from WASP-33 

4.1. Brightness and equilibrium temperatures 


From the observed planet-star flux ratio at a given wavelength 
the exoplanet brightness temperature Tb at that wavelength can 
be derived. For this calculation we represent WASP-33’s emis¬ 


sion by a spectrum from the PHOENIX library (Hauschildt 


|& Baron] 1999| [Witte et al.| |2009| l with stellar parameters 
Teff = 7400 K, [Fe/H] = 0.1, and log g = 4.5, closely matching 
WASP-33’s values given by Collier Cameron et al. ( |2010| l. 
The planet contribution is represented by a black body and 
the measured eclipse depth represents the planet-to-star flux ra¬ 


tio Fp /F,s(Y). integrated within the F-band (Charbonneau et al. 
2005 1 . We computed between our best-fit eclipse depth and 


Fp/Fs(Y), varying Tb between 2500 and 4500 K each 1 K. By 
means of x^ minimization we found the brightness temperature 
of the exoplanet to be Tb(Y) = 3398 + 302. Its error was deter¬ 
mined from a Ax^ contour equal to 1 . 

To date, six measurements of WASP-33 b’s eclipses along 
with their corresponding have been obtained (Table [^. 
Assuming that the brightness temperatures can represent the 
equilibrium temperature Teq of WASP-33 b (see e.g., de Mooijj 
jet al.||2013] Section 4.1), we re-determined it following the ap- 
proach explained in the previous paragraph but using all bib¬ 
liographic data instead. We found the best-fit equilibrium tem¬ 
perature equal to Tgq = 3358 + 165 K (Figure p]i. Our determi¬ 
nation of the equilibrium temperature is ~ 10 frK higher than, 
but co nsistent with, the value determined by jde Mooij et al.j 
|2013 I, Teq = 3298)^g® K. The corresponding Yjgd ^ d.o.f. is 
1.2. We calculated values between the observations and the 
two cases presented in [Deming et al.| ( [2012| l (see their Figure 9 
and our Figure]^, i.e. a carbon-rich non-inverted model (green, 
X^^^ - 1.95, 3 d.o.f.), and a solar composition model with an in¬ 
verted temperature structure (pink,;if^^j = 5.71, 3 d.o.f.). We find 
that our simple blackbody curve seems to provide comparable 
results to the non-inverted scenario, and could also easily fit the 
inverted scenario if the emission from the TiO band centered at 
0.9 pm is stronger than what current models predict. We notice 
that the measurement at 2 microns differs significantly with re¬ 
spect to the inverted (pink) model. This discrepancy has been 
solved by [Haynes et al.[([2015[) (see their Section 5.1 and their 
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Figure 5). We believe that these data can not be uniquely repre¬ 
sented by a given model without the need of further speculation, 
and find our simple blackbody approach as sufficient. A more 
sophisticated analysis using different models with different in¬ 
verted temperature structures is beyond our scope. 

Special notice has to be given to the wavelength coverage of 
the data. The lower panel of Figure|^shows that all planet-to-star 
flux ratio measurements lie within the Rayleigh-Jeans region of 
the planets black body, making any temperature measurement 
prone to errors, since the raising part of the planetary emission 
curve is not being taken into account. A flux measurement at 
A ~0.7-0.9 yum would significantly improve our knowledge on 
WASP-33 b’s temperature, because the shape of the planetary 
SED would be better constrained and the determination of the 
planetary temperature wouldn’t be determine only by the slope 
of the SED. However, secondary eclipse observations in the op¬ 
tical are extremely challenging, aggravated in this particular case 
by WASP-33’s intrinsic variability. 


Table 5: Observed eclipse depths (EDs) of WASP-33 b in ppt 
and the derived brightness temperatures (Tb). CW corresponds 
to the central wavelength of the filters. 


Authors 

EDs 

Tb 

CW 


(ppt) 

(K) 

(/tm) 


Smith et al. 

2011 



1.09 + 0.30 

3490 + 140 

0.91 

Deming et a 


( 

m 

I'l 

1 

2.70 + 0.40 

3415 ± 130 

2.14 

Deming et al. 

( 

2017 

1 

2.60 + 0.50 

2740 ± 225 

3.55 

Deming et al. 

( 

2012 

1 

4.10 + 0.20 

3290 ± 100 

4.50 

de Mooij et 

a 


( 

wr 


2.44 + 0.25 

3270+115 

2.14 

Ihis work 






1.03 + 0.34 

3398 + 302 

1.05 


4.2. Albedo and re-distribution factor 

Due to their proximity to the host star, the exoplanet atmospheres 
of short period exoplanets are strongly dominated by the stel¬ 
lar incident flux. An indicator for the fraction of light that gets 
reflected back into space immediately is the geometric albedo 
for which the majority of measurements yielded very low values 
( |Heng & Demory||2013[ ). How well the absorbed flux is redis¬ 
tributed from the dayside to the nightside is still a matter of de¬ 


bate (|Pema et al. 

2012||Heng & Showman|20141i.|L6pez-Morales| 

|& Seager (2007 

1 presented an analytic expression for the exo- 


planet equilibrium temperature as a function of the Bond albedo 
Ab and the reradiation factor /, which describes how the stellar 
radiation absorbed by an exoplanet is redistributed in its atmo¬ 
sphere. Using their Eq. 1, we compute the equilibrium tempera¬ 
ture of WASP-33 b for a zero-albedo and reradiation factors of 
2/3 (instant reradiation) and 1/4 (instant redistribution) obtaining 
Teq = 3499 K and T^q - 3257 K. These are shown in the bottom 
panel of Eigurej^as orange diamonds and cyan triangles, respec¬ 
tively. 

Instead of choosing two boundary conditions for Ab and /, to 
constrain WASP-33 b’s values of bond albedo and redistribution 
efficiency we produced an equilibrium temperature map. Each 
point of the map was calculated in the same fashion as before, 
but defining a 300x300 grid in the Ab - f space. The map is 
displayed in Eigure|^ where the cyan area defines WASP-33 b’s 
equilibrium temperature range. The shaded effect was produced 
considering the stellar (Tefj, Rs) and orbital (a/Rs) parameters 
with their uncertainties. 



0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 


Wavelength (^m) 


Pig. 4: Top: Eclipse depth in ppt as a function of wave¬ 
length. The blue data point indicates our brightness tempera¬ 
ture measurement, while the red circles show literature values. 
The gray area shows 1-cr contour of WASP-33 b’s equilibrium 


temperature. Green continuous line corresponds to Deming et al. 
( |2012| l’s carbon-rich non-inverted model, while pink continuous 
line indicates a solar composition model with an inverted tem¬ 
perature structure. Squares mark the planetary fluxes integrated 
over the bandpasses. Bottom: The same data after multiplication 
by the PHOENIX stellar model revealing the Spectral Energy 
Distribution (SED) for WASP-33 b. Overplotted are the expected 
SEDs for a zero-albedo instantly re-radiation (orange diamonds, 
f = 2/3) and instantly redistribution (cyan triangles, f = 1/4) day- 
sides. 


Producing a statistical analysis of the bond albedo and re¬ 
distribution efficiency over 24 transiting exoplanets with at least 
one published secondary eclipse, |Cowan & Agol| ( |201 1 1 found 
that all the planets in the sample were consistent with low 
bond albedos (< 10%), and planets hotter than 2400 K pre¬ 
sented low redistribution efficiency (i.e., large / values). This 
is in good agreement with the cyan contours displayed in Eig- 
ure Thus, WASP-33 b seems to have similar properties as 
atmospheres of other close-in gas giants orbiting cooler stars 
(Ab <0.38, 0.41<f<2/3). 

A nonzero geometric albedo of the planet in the Y band 
would result in a contribution of reflected light to the eclipse 
depth. Eor example, the potential impact of a rather high geomet¬ 
ric albedo of 0.5 would be ~350 ppt on the eclipse depth, which 
is on the order of our 1-cr uncertainty. Such contribution of re¬ 
flected light would not significantly affect our results of Sect. 
4.1. However, any nonzero albedo would lower the estimated 
exoplanet equilibrium temperature. 
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Fig. 5: Equilibrium temperature map as a function of albedo 
(Ab) and re-distribution factor /. Cyan color shows the contours 
relative to WASP-33 b’s equilibrium temperature within 1-cr er¬ 
ror (top-left corner) for both planetary and stellar parameters, 
f values larger than 2/3 are unphysical. Therefore, the white- 
shaded area of the plot has not to be considered. 


4.3. Evidence for temperature inversion? 


The observed diversity in the characteristics of hot Jupiter atmo¬ 
spheres has been a main subject of investigation. For instance, 
Knutson et al.| ( T010| l suggested a direct relation between the 
activity of the host star and the presence of a temperature inver¬ 
sion. While active stars should host exoplanets without temper¬ 
ature inversions, the opposite should occur around quiet stars. 
To investigate WASP-33’s activity levels we analyzed calcium 
lines in high-resolution spectra acquired on December 28, 2012, 
using the High Resolution Echelle Spectrometer mounted at the 
Keck 1 telescope (Program ID N116Hr, PI Howard). WASP-33 is 
X-ray dark and generally inactive, as seen in a lack of emission 
in the Ca II H and K lines (Figure [^. Comparing the observa¬ 
tions to [Knutson et al.| ( 2010[ )’s suggestions, WASP-33 b should 
present a temperature inversion. Also, planets with inversions 
are expected to have greater SPITZER/IRAC 4.5 pm than 3.6 
pm flux, resulting from strong CO and CO 2 absorption in the 
4.5 pm bandpass (see e.g., [Burrows et al.||2007{ [Madhusudhan[ 
" ' 2010[ ). This has been observed in WASP-33. However, 


& Seager 


both results seem t o be in disagreement with the observed high 
re-radiation factor (Fortney et al.[2008i. Although analysis 


would favor a non-inversion, the here presented disagreements 
caution ourselves to make any assumption about the atmospheric 
structure of WASP-33 b. A phase curve would help to constrain 
the albedo and the re-radiation factors. However, stellar pulsa¬ 
tions might make the task more challenging. 


5. Conclusions 

In this work we present simultaneous observations in the V and Y 
bands of WASP-33 during a secondary eclipse. We make use of 
the information about the stellar pulsations on the V light curve 
to “clean” the NIR data, minimizing the impact of stellar pul¬ 
sations and correlated noise in the determination of the eclipse 
depth. The measured eclipse depth is AF = 1.03 + 0.34 ppt, 
which corresponds to a brightness temperature of 3398 + 302 K. 
We represent the spectral energy distribution of WASP-33 b by 
a black body curve and fit it to current data. We find an equi¬ 
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Fig. 6: WASP-33’s high resolution spectra around the Ca II H&K 
lines observed with Keck/HIRES. 


librium temperature equal to 3358 + 165 K, slightly hotter but 
still consistent with previous measurements. Comparing val¬ 
ues between our approach and more sophisticated measurements 
we conclude that these observations can not provide a detailed 
description of the temperature structure of WASP-33 b’s atmo¬ 
sphere. In agreement with previous publications and hypothesis 
made by Cowan & Agol| ( [201 Ij l, we find that WASP-33 b has 
a very low efficiency to circulate stellar radiation into the exo¬ 
planet night side. 
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